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Ideal magnetohydrodynamics (IMHD) is strongly constrained by an infinite number of 
microscopic constraints expressing mass, entropy and magnetic flux conservation in each 
inhnitesimal fluid element, the latter preventing magnetic reconnection. By contrast, in 
the Taylor relaxation model for formation of macroscopically self-organized plasma equi¬ 
librium states, all these constraints are relaxed save for global magnetic fluxes and helic- 
ity. A Lagrangian variational principle is presented that leads to a new, fully dynamical, 
relaxed magnetohydrodynamics (RxMHD), such that all static solutions are Taylor states 
but also allows flow. By postulating that some long-lived macroscopic current sheets can 
act as barriers to relaxation, separating the plasma into multiple relaxation regions, a 
further generalization, multi-region relaxed magnetohydrodynamics (MRxMHD) is devel¬ 
oped. 

PACS codes: Authors should not enter BAGS codes directly on the manuscript, as these 
must be chosen during the online submission process and will then be added during the 
typesetting process (see http://www.aip.org/pacs/ for the full list of PACS codes) 


1. Introduction 

The coarse-grained, fluid-like dynamical behaviour of highly conducting, magnetized 
plasmas in the laboratory, in stars such as the sun, and in space, can often be described 
by variants of magnetohydrodynamics (MHD). 

The one-fluid, non-dissipative model most commonly used, ideal magnetohydrodynam¬ 
ics (IMHD), appears deceptively simple but is strongly constrained by an infinite num¬ 
ber of microscopic constraints expressing the detailed conservation of mass, entropy and 
magnetic flux “frozen” into each infinitesimal fluid element. Physically, these are very 
restrictive constraints, e.g. the “entropy freezing” constraint prevents dynamical temper¬ 
ature equilibration along magnetic field lines, and the “flux-freezing” constraint prevents 
changes in the topology of magnetic field lines, thus preventing magnetic reconnection, 
island formation, or formation of chaotic lines. These constraints also give rise to math¬ 
ematical problems due to a tendency for singularities to form when systems are per- 


turbed away from simple geometries with a continuous symmetry. 
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Figure 1. Some possible relaxation regions, discussed in Sec. |2.2| The innermost region, fio, 
is a simple genus-1 toroid bounded by a torus dflo, formed from Fj^j^ and (blue and 


orange), the in ward faces of the inner separatrix current sheets, Waelbroeck (19891; [Wang fc 


Bhattacharjee (1995 l, of an m = 2, n = 1 magnetic island chain. The helical interiors of the two 
islands of the chain, flia and Slij, are also genus-1 toroids, bounded by the tori obtained by 
joining with F^j, S'lid with rj"j, where j a,nd Fjj, j ®re the inward faces of the 

outer separatrix current sheets (blue and orange). (Colour online.) 


to formulate an alternative self-consistent, non-dissipative single-fluid model for toroidal 
plasmas that is simpler than IMHD yet is less physically restrictive and better posed 
mathematically in general geometries. 

The most powerful and general way to formulate a non-dissipative field theory, see 
e.g. p. 53ff of Goldstein (1980), is to postulate a Lagrangian density £ and to derive 


the dynamical equations for all fields from the action, S = f dtf^ £ d^x, by appealing 
to Hamilton’s Principle. That is, by requiring that its first variation, SS, vanish for all 
variations of the independent fields in the system region G. The equations for these 
fields are the resulting Euler-Lagrange equations, which, as they are all derived from 
the one scalar functional S, are automatically self-consistent. Furthermore, conservation 
equations can be derived very generally by applying Noether’s Theorem [see e.g. p. 555ff 
of Goldstein (1980), also K. Charidakos et al. (2014) and references therein], based on 
the continuous symmetries of the system. Holonomic constraints can be handled by 
expressing variations of dependent fields in terms of those of the independent fields, 
and non-holonomic constraints by augmenting the Lagrangian density using Lagrange 
multipliers. Our modified MHD is based on the same Lagrangian as IMHD, but uses a 
much-reduced set of constraints, a small subset of those implicit in IMHD. 

In IMHD and our modifications of it, on the boundary 9G of the overall plasma region 
(and on current sheets separating plasma subregions), the magnetic field B is constrained 
everywhere to be a tangent vector 


n ■ B = 0 


(I.l) 


where n is the unit normal at each point on dil. Physically, this corresponds to the 
assumption of confinement within a perfectly conducting wall [but not necessarily a rigid 
wall if one wishes to model, for example, the response to an externally imposed pertur¬ 


bation by switching on boundary ripple, Hahm & Kulsrud (1985); Dewar et al. (2013); 
Comisso et al. (2015 In this paper we assume for simplicity that the wall has no 
gaps, so that the wall completely shields the plasma from penetration of externally gen¬ 
erated magnetic fluxes, but this restriction is not essential for a Lagrangian formulation 


to be possible, Dewar (1978); Hosking & Dewar (2015), and would need to be lifted if 


one wished to consider Ohmic current drive or helicity injection. 

In the following we will partition into interacting subregions £li, on the boundaries of 


which the tangential-B constraint eq. (1.1) is also enforced, but tangential discontinuities 
due to current sheets on these interfaces are allowed. (A vacuum region Gv in which no 
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IMHD invariants other than total magnetic fluxes are assumed, can also be included 
between the wall and a plasma-vacuum interface.) 

We also make the topological assumption that the plasma regions 17^ are toroids, so 
their boundaries dQ are tori, where by (generalized) tori we mean 2-dimensional surfaces 
without boundary, conceivably multi-handled, and by toroids we mean 3-dimensional 
volumes bounded by a single torus (which we term a simple toroid), or by an inner and 
an outer torus (an annular toroid). (See Fig. [^) 

The single-fluid, non-dissipative MHD equations are encapsulated in the Lagrangian 


density, Newcomb (1962); Dewar (1970), 


rMHD ^ ^^^2 


7-1 


2lio 


( 1 . 2 ) 


where p is the mass density, p is the pressure, and po is tbe permeability of free space. 
The plasma is treated thermodynamically as an ideal gas with an isentropic (adiabatic) 
equation of state, p/p'^ = const (in IMHD this being applied microscopically in each fluid 
element). 


The ideal (IMHD) equations of motion follow, Newcomb (1962); [Dewar (1970), from 


eq. (1.2) by applying Hamilton’s principle of stationary action, treating trial displace¬ 


ments of fluid elements from their physical positions as an arbitrarily variable vector 
field, to which the variations of density p, and pressure p (or entropy density) and mag¬ 
netic held B, are holonomically constrained microscopically. A (non-canonical) Hamil¬ 
tonian formulation of IMHD, in which constraints appear as degeneracies of Poisson 


brackets (some of which can be integrated as Casimir invariants), is also possible, Mor- 
( I998|); [Yoshida fc Dewar (2012), but the Lagrangian approach provides a more 


rison 


convenient starting point for hnding a modihed magnetohydrodynamics. A Hamiltonian 
formulation could be derived from our new Lagrangian formulation, but this is not pur¬ 
sued in the present paper. 


In strong contrast to IMHD, Taylor’s relaxed equilibrium model, Taylor (1974), re¬ 


laxes all the IMHD constraints save for conservation of toroidal magnetic hux and global 
magnetic helicity (which are IMHD invariants), leading to a very low-energy “relaxed” 
equilibrium state. Such a “Taylor state” is a special static solution of the IMHD equa¬ 
tions, but is dynamically inaccessible via IMHD from arbitrary initial states because of 
IMHD’s inhnity of extra constraints. To elevate Taylor’s static relaxed equilbrium theory 


to a relaxed magnetohydro-(RxMHD), we use the Lagrangian eq. (1.2) and 
the same holonomic density constraint as for IMHD, Newcomb (1962); Dewar (1970), 


but treat the pressure and magnetic held as independently variable helds subject only 
to conservation of global hux(es), entropy and magnetic helicity within D, and the holo¬ 
nomic tangential-B constraint eq. (1.1). This makes the Taylor equilbrium state always 


dynamicallly accessible within RxMHD. 

Physically, the Taylor model is designed to predict the hnal macroscop*([f]self-organized 
state to which a highly conducting plasma will evolve, provided it has a stochastic mech¬ 
anism for breaking the “freezing in” of magnetic held by the plasma, [Rusbridg^ (1991); 
Qin et al. (2012), and for generating magnetic held by small “dynamo” hows. Since Tay¬ 


lor’s pioneering work, Taylor (1974), which invoked the conservation of magnetic helicity 


I By “macroscopic” we mean a coarse-grained description in which the small-scale reconnec¬ 
tion processes are spatially and temporally nnresolved, and also a description that is only correct 
in the limit that the non-dimensional conductivity parameter (Lundquist number) increases to¬ 
ward infinity, the dynamo flows being assumed to be higher order in inverse Lundquist number 
so that the Taylor state is a static, force-free MHD equilibrium. 
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in the relaxation of toroidal discharges, an extensive literature has arisen on this topic, 


much of it reviewed by Taylor (19861. 


Taylor’s original theory, Taylor (1974), applied globally throughout the plasma with 


very few adjustable parameters, and was remarkable for its success in modelling toroidal 
held reversal and helical bifurcation in the highly turbulent reversed-held pinch, Zeta. 
However, this simplicity restricts its ability to model better-conhned axisymmetric plas- 
Bhattacharjee & Dewar (1982), such as tokamaks or more modern reversed-held 


mas, 


pinches, whose modelling needs more constraints to increase hexibility in matching ob¬ 
served prohles. It is even less adequate for modelling non-axisymmetric systems, such as 


tokamaks [see e.g. Figs. 7 and 8 of Hudson et al. (2012)] or reversed-held pinches, Dennis 
et al. (20136), with weakly broken symmetry; or stellarators designed from the outset to 


be non-axisymmetric. 


As will be explained further in Sec. M we are led to generalize Taylor’s theory by 
replacing the smooth constraints of Bhattacharjee & Dewar (1982) with singular con¬ 
straints —macroscopic current sheets, Fjwhich partition the plasma into multiple re¬ 


laxation regions [These current sheets may be thought of as thin, hexible sheets of 
ideal plasma, within which all the IMHD invariants apply, corresponding to the “singular 
Casimir elements” of Yoshida & Dewar (2012).] This leads to a further dynamical gener¬ 


alization, multi-region relaxed MHD (MRxMHD), which we anticipate will have a number 
of applications in modelling toroidal confinement devices in which three-dimensional ge¬ 
ometry effects are important. Current-sheet and relaxation theory applications are also 


ubiquitous in astrophysical contexts, Parker (1994), though one has to deal in these 


applications with the constraint of line-tying, not covered in this paper. 


A new numerical approach to calculating plasma equilibria, using a static version of 
MRxMHD, has been implemented in a code, SPEC, Hudson et al. (2012), particularly 


useful in non-axisymmetric toroidal fusion confinement systems when most flux surfaces 
are destroyed by held-line chaos but also capable of reproducing ideal-MHD calcula¬ 
tions in systems with a continuous symmetry by using a large number of nested annular 


toroidal subregions, Dennis et al. (2013a). This static formulation has also recently been 
used to explore MHD singularities at resonant magnetic surfaces, Loizu et al. ( 2015a|6 ). 
Interestingly, our MRxMHD formulation has some striking similarities to an early com¬ 
putational “water bag” approach proposed, but apparently not developed further, by 


Potter (1976). 


This paper attempts to construct a general formal framework for MRxMHD from 
first principles in a pedagogic manner, citing related historical and recent work where 
possible, indicating the scope of MRxMHD, and setting the stage for further development 
and application. Taylor relaxation and helicity conservation are reviewed in Sec. 


2.1 


with generalizations discussed in Sec. |2.2| In Secs. |3.1| and |3.2| we review standard fluid 
Lagrangian variational results, using a notational framework that is formally precise for 
use in further work. The Lagrangian framework for MRxMHD is developed in Sec. |3.3[ 
and MRxMHD dynamics is derived from Hamilton’s Principle (of stationary action) in 


Secs. 3.4 and 3.5 with constraints corresponding to our chosen subset of IMHD invariants. 
Possible further developments are suggested in the Conclusion. Appendix reviews 
the simple thermodynamics used in MRxMHD and Appendix reviews the boundary 
condition for the vector potential. 
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2. Plasma relaxation 

2.1. Helicity conservation and Taylor relaxation 
The conservation of magnetic fluxes threading cuts through (or 12^ in MRxMHD) 
that leave it topologically connected follows simply from V • S = 0 and the boundary 
condition eq. (O everywhere on dTl (because of the no-gaps assumption mentioned 
in the Introduction). Thus, while these fluxes are IMHD invariants, they must also be 
invariants in any physical model — flux conservation is said to be completely robust. 

As it is conserved even under reconnection [in the limit as resistivity approaches zero, 
see e.g. eq. (20) ff. of Jensen & Chu (1984|] and has a robust topological interpretation 
[see e.g. Berger ( 1999| ) for a heuristic review or Arnold fc Khesin] (1998) for a more 
mathematical treatment] the most robust of the remaining IMHD invariants is widely 


accepted to be the magnetic helicity 2]1 qKq, where, Bhattacharjee & Dewar (1982), we 
define the invariant Kq as 

f A-B , , 

Kn= -dH, (2.1) 

with A a vector potential giving B — 'VxA and dV = d^a; the volume element. Because 
K has one less gradient of A than the magnetic energy. 


W^ = 


2po 


-dB 


( 2 . 2 ) 


it can also be argued, see e.g. Sec. I.C of Taylor (1986), that, in a weakly resistive plasma 
with small-scale turbulent fluctuations, K decays slower with time than . 

Other general ideal invariants can be related, by Noether’s theorems, to some symme¬ 
tries in appropriate parameterizations of field variables in the IMHD action (for example, 
the cross helicity f^v • BdV pertains to a relabelling symmetry in the Lagrangian rep¬ 
resentation of the fields), Salmon (1988); Padhye & Morrison (1996 |Webb fc Zank| 
(2007); Webb et al. ( 2014a\b ); Araki (2015). However, in the spirit of Taylor relaxation 
we choose the minimal set required to obtain a non-trivial solution and thus keep only 
Kfi as the only non-holonomic constraint involving v or B. 

It is readily shown, using the tangential-B condition eq. (1.1), that Kq is invariant 


under gauge transformations A i-^ A -|- Vx as long as x is single-valued (implying con¬ 
servation of line integrals /g^A • dl around loops on the boundary, which, by Stokes’ 
theorem, is equivalent to the above-mentioned conservation of magnetic fluxes). As tan¬ 
gential B is to be a holonomic constraint rather than a natural boundary condition, we 
do not treat A on the boundary as freely variable and can constrain y. Thus we do not 
need to use either the Bevir-Gray [subtraction of products of toroidal and poloidal loop 
integrals, Bevir fc Gray] ( |l9^ )] or relative helicity (subtraction of vacuum-field helicity, 
Jensen & Ghu (1984)) modifications of the helicity, the latter fact also implying there 
is no physical necessity to decompose the magnetic field into a vacuum (harmonic) and 
a plasma-current-generated component [though it may still be useful conceptually and 
mathematically, Yoshida & Giga (1990); Yoshida & Dewar (2012)]. 

The Woltjer-Taylor variational principle [originally proposed, though with less physical 
motivation, by Woltjer (1958)] is that the final relaxed state is that which minimizes the 
magnetic energy, eq. (2.2) [the negative of which occurs in eq. (1.2)], under the magnetic 
helicity constraint, implemented by minimizing Wq — under variations of the mag¬ 
netic vector potential A, fi being a Lagrange multiplieijt] The resulting Euler-Lagrange 


I Such Beltrami constants p have dimensions of inverse length and are not to be confused 
with the vacuum permeability constant po used in SI units. 
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equation is the Beltrami equation, 

VxB = fiB. (2.3) 

This describes a force-free field, i.e. one with the current j = V X B/po parallel to B, 
implying Vp = 0 globally in an equilibrium plasma. Thus this single-region relaxation 
principle describes only plasmas with no thermal confinement. 


2.2. Generalization of Taylor Relaxation 

As the Woltjer-Taylor variational principle is not explicitly based on knowledge of the 
detailed sub-macroscopic physics leading to relaxation and self-organization, its applica¬ 
bility to modelling a given system can only be justified empirically|f] We have already 
remarked in the Introduction that in fact Taylor relaxation theory in its original form 
is too simple to apply to modern fusion devices, but that it can be extended naturally 
by supplementing magnetic helicity with further global ideal invariants, thus preserving 
much of its simplicity but increasing its flexibility in applications to modelling fusion 
plasmas. 

The fundamental basis of our generalized MHD relaxation principles is the requirement 
that the states they describe be a subset of the states allowed within ideal MHD. This is 


ensured, Bhattacharjee & Dewar (1982); Dewar et al. (2008), by using only constraints 


from a subset of those implied by IMHD. We take this formal criterion as the paramount 
principle for constructing consistent modifications of IMHD, regarding considerations of 
possible sub-scale physics that might lead to breaking of some ideal invariants and not 
others only as an heuristic guide in choosing an appropriate subset of IMHD constraints. 
For instance, in choosing relaxation subregions in which to apply the Woltjer-Taylor 
variational principle we do not necessarily assume the magnetic field is wholly or partially 


chaotic, though Beltrami solutions can accommodate such cases, Dombre et al. (1986). 


Justification for the choice of constraints must ultimately be empirical, by comparison 
either with experiment or ab initio simulations. 

A generalization of the Taylor relaxation idea by increasing the number of constraints 


was proposed by Bhattacharjee & Dewar (19821, but the smooth IMHD invariants chosen 


then are not well-defined in a nonintegrable magnetic field with islands and chaotic 
regions. More recently, Hudson et al. (2007); Dewar et al. (2008), generalizations of 


Taylor relaxation theory were proposed, based on the assumption that Taylor-relaxed 
plasma can coexist with current sheets that act as transport barriers partitioning Q into 
multiple regions invariant under field-line flow, Hudson et al. (2012). To describe 


this approach we have introduced the terminologies RxMHD when fl is not partitioned, 
and MRxMHD when it is, the D (for “dynamics”) being justified below. The J-function 
currents in MRxMHD are compatible with IMHD so they may be regarded as singular 
alternatives to the smooth IMHD constraints of Bhattacharjee & Dewar ( |1982 ). As in 
[Bhattacharjee &: Dewar (19821 we also introduce entropy constraints to allow a nontrivial 


pressure profile and retain the non-singular magnetic helicity invariant(s), eq. (2.1), but 
separately conserved in each MRxMHD subregion. 

Our development of MRxMHD is implicitly based on the multiple-timescale scenario 
sketched below (a conceptual framework motivating the formal development—whether 

I This indeed is also true of ideal MHD, which is typically applied in fusion physics well 


beyond the validity of the 
e.g. Sec. II.H of Freidberg 
field lines are not short in 
parallel and perpendicular 

approximations required for using it to describe fusion plasmas [see 
(1982l]. In particular, particle mean free paths parallel to magnetic 
iigh-temperature plasmas, so modifications of IMHD that distinguish 
physics have long been sought in order to extend its applicabilty [e.g.. 

the collisionless MHD of 

Freidberg 

(1987 

H- 
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there are situations where it approximates physical reality remains to be explored in 
further work): 

1. A fast relaxation timescale during which all but a finite number of IMHD constraints 
are broken through thermal diffusion and micro-reconnection events (associated with 
unspecified mechanisms like micro-tearing turbulence, high-order resonant structure^ 
and field-line chaos). On this timescale the system self-organizes into multiple Taylor 
states in disjoint subregions with non-disjoint boundaries dfli (geometrically fixed 
on this timescale due to plasma inertia) supporting current sheets on their common 


interfaces Tij [recent simulations by 

Smiet et al. 

(^ 

>015 

) give some support for this 

scenario]. The tangential-B boundary condition, eq. 

(1.1 

) is satisfied on both sides of 


these interfaces, but in general B suffers a tangential discontinuity across them. [To 
represent discontinuities across an interface T^ j, we distinguis h its inward and outwa rd 
faces Tfj by the superscripts — and -f, respectively, following McGann et al. (2010).] 


2. An intermediate dynamical timescale (the timescale treated in this paper) during 
which the plasma, including a number of embedded current sheets separating sub- 
regions within which magnetic helicity, magnetic fluxes, mass, and entropy are con¬ 
served, evolves adiabatically with respect to the relaxation timescale as the geometric 
shapes of the boundaries dfli evolve dynamically from their initial conditions, and 
possibly in response to external forcing from the “switching on” of boundary ripple. 


Hahm & Kulsrud 

(1985 

I; 

Dewar et al. (2013 

); Comisso et al. 

(2015a|l))). Low-order 

resonant structures within the plasma that are excited. 

Boozer & Pomphrey 

(2010); 


White (2013), by geometric change resist the formation of magnetic islands by devel¬ 


oping shielding current sheets. 


3. A long reconnection timescale on which plasma and magnetic flux leaks and mixes 
between sub-regions through weak spots in the current sheets Ti^-, violating the mass 
and flux isolation of the sub-regions assumed in MRxMHD and also violating entropy 
conservation. New subregions may form, changing the topological structure of the sys¬ 
tem. Phenomena on this timescale are not treated in this paper. 


To illustrate these concepts, a multi-region case of interest is shown schematically in 
Figs.[^and[^ an island chain {flia, twisting, with linking number 2, Berger (1999), 
around a simple toroid Dq, both contained within an annular toroid 122- This illustrates 
the flexibility of the MRxMHD partition—it is not limited to simply nested tori as it can 
include island or plasmoid chains. (This requires a more complicated labelling system 
for the relaxation regions and current sheets, but if we limit ourselves to primary island 
chains we can still use a sequential numbering system and indicate the component islands 
and separatrices of the chain using the lettering scheme illustrated here.) 

In single-region RxMHD the topology of determines the number, of independent 


Berger (1999 


t We use the terminology “resonant structure” to denote a family of closed field lines, whose 
order is the number of toroidal rotations they make before they close. A resonant toroidal flux 
surface is a special case of s uch a structur e, but transient resonant structures may form due to 
Sweet-Parker reconnection, Parker (19941, initiated at initially isolated hyperbolic closed field 
lines s u ch as the “X points” of mal efic islands. While numerical evidence, [Longcope fc Strauss] 
(19931; Cordoba & Marliani (20001, that strictly d-function current sheets can form in finite time 
is not conclusive, it is a reasonabie postulate in our coarse-grained, long-reconnection-timescale 
MRxMHD model. 

j: The topologically genus (or hrst Betti number) v is the number of cuts required to make 
D simply connected. E.g. in a simple toroid (cf. Do in Figs. and [^, i/ = 1, as a toroidal cut 
leaves it simply connected. 
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Z f = 0 Z f = ;r/2 



Figure 2. Four sections of the regions shown in Fig. at toroidal angles ^ = 0, 7r/2, tt and 
37r/2, showing helical rotation of the islands ilia and idu about the central region flo- The outer 
region, 0 , 2 , is an annular genus-2 toroid bounded by F^^ ^ s-nd the plasma-vacuum 

interface (green torus). 


fluxes \ ^v, that must be specified for uniqueness of a harmonic (vacuum) held 

Bh solution such that VxSH = 0inn,n-BH = 0oni9n. For not an eigenvalu^of 
th e Beltrami equ ation, eq. (2.3), with homogenous boundary conditions [see e.g. Sec. IV 
of Taylor (1986)], specifying the u fluxes also specifies the Beltrami field uniquely, 


Yoshida & Giga (1990). Similarly, in MRxMHD we need to determine the genus of 
each relaxation sub-region and specify its fluxes <!>[, which are invariant under relaxation. 

In the case of an annular toroid Oi (e.g. O 2 in Figs. [^and[^, dOi consists of two 
disjoint tori (e.g. Fj*)j UFj*), and F^). A standard single-handled torus (i.e. with one hole) 
can be covered by a single coordinate chart, typically using a poloidal angle 9 and a 
toroidal angle C. Assuming its boundaries to be two such standard tori, the genus of an 
annular toroid is u = 2 , as both toroidal and poloidal cuts are required to make it simply 
connected. 

Higher-genus cases might arise in toroidal confinement when treating doublet/multi¬ 
pinch, Taylor (1986 1 , or bundle divertor, Stott et al. (1977), configurations. In the bundle 
divertor case, the plasma-vacuum boundary dOl is a two handled torus, which cannot be 
described by a single toroidal-poloidal coordinate system but must instead be partitioned 
into two separate patches (coordinate charts). We shall not consider such exotic cases 
further in this paper. 

The conservation of the <I>[ implies boundary constraints on the vector potential A: 
From Stokes’ theorem 


B-ndS = j^A-dl (2.4) 

where n is the unit normal at a point on the Zth topologically distinct surface of section 
(t[ cutting Oli and 7 * = da{ S dOli is a loop around the boundary of cr[, with direction 
with respect to that of n given by the right-hand rule. 

As the loops 7 j- lie on the boundary dOti, which is composed of the current-sheet 
interfaces Fi ^ where there are 5-function currents causing tangential discontinuities in 
B, one might think the values of the loop integrals would depend on whether the loops 
traverse the inner or outer faces of the F^j . However, it does not matter which faces are 
used as B remains finite within the current sheet—hence, being of infinitesimal width, a 



^ In the present context this means gi below the lowest Beltrami eigenvalue in each Oi, else 
the plasma would be unstable to local tearing instability, in which case Oi should be partitioned 


further to raise the minimum eigenvalue. 


Dennis et al. 


(20136). 
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current sheet contains only infinitesimal flux. However, this continuity of its loop integrals 
does not necessarily mean A itself is continuous across current sheets, as a discontinuous 
gauge term Vy does not affect the loop integrals (provided x is single valued). This 
freedom allows a coordinate-dependent gauge to be used in each relaxation region, as in 
the SPEC code, Hudson et al. (20121. 

It is important to recognize that the fluxes <!>( depend only on magnetic fields within 
Qi, so their individual conservation constrains line integrals of A only around loops 
that enclose the plasma within However, taking into account the conservation of all 
the leaves only the toroidal line integral fg^A • dl on the plasma-vacuum interface 
unconstrained. This represents the external poloidal flux threading the hole in the torus 
While this flux is arbitrary as far as the physics of the plasma within fl is concerned, 
it is still conserved because we are assuming the wall acts as a superconducting shell 
which traps the external poloidal flux threading it. 


3. Lagrangian formulation via Hamilton’s Principle 
3.1. Lagrangian and Eulerian fluid kinematics 

Central to the Lagrangian approach to fluid mechanics is the concept of fluid elements, 
whose motions with respect to time t through a 3-dimensional Cartesian frame (points 
in which we designate by the vector x = xe^ + ye.y + zeflj, are described by a family of 
feasible trajectories (pathlines respecting the constraints) a:* = r‘(a;o), labelled by Xq, 
the initial positions of fluid elements at an arbitrary time t = to. The fluid elements 
are advected by the Eulerian velocity vector field v{x,t) to their positions at arbitrary 
time t through the time evolution function r* (x) defined as the integral of the following 
equation and initial condition 


dr*(a;) 

dt 


= v[r 


{x),t) 




{x) = X 


(3.1) 


for all X in the domain of interest. (Note that this makes r* implicitly a function of toi 
which we can make explicit when needed using the notation r‘(a:|to), constant parameters 
such as to being listed after the vertical bar |.) 

Where an argument (other than t) is not specified, we treat r* as a map —>■ (i.e. 
a 3-vector function of 3-vectors), but often it is necessary to recognize that it is also a 
functional of v{x,t), which will be indicated explicitly when required using the notation 
r^[v]{x). This provides a very flexible notation that can be adapted for generating other 
maps. 

Suppose V is also a function of some time-independent parameter [say s, denoted 
by t;(a;,t|s)] then r* will also be a function of s, denoted r^{x\s,to). Suppose further 
that X lies on an arbitrary curve x = f{s) at time to and denote the resulting family 
of trajectories by R*{s) = r‘(/(s)|s, to)- Differentiating both sides of the equation of 


motion in eq. (3.11 with respect to s, we have 


d dR* dRf „ ^ 

• Vv{R*,t) 

dt ds ds 


dsv{R\t) 


(3.2) 


where ds means the partial derivative with respect to s and Vv denotes 'Vxv{x,t). 

We now use the special case of the above result where v is independent of s to build 
up some useful differential-geometric evolution resul ts. F irst, denoting an infinitesimal 
line element advected by the fluid by dl* = dR*, eq. (3.2) immediately gives 


^dt* = dl* • Vv . 
dt 


(3.3) 
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When acting on all points in a region il, the evolution function r* defines the La- 


grangian map, mapping an initial region fig onto its image fl*. From eq. (3.3) applied to 


the sides of an infinitesimal rhomboid within the advection equation for infinitesimal 
volumes dF* is found to be 

d 


dt 


dW = (V • t;)dW 


(3.4) 


which is equivalent to the evolution equation, (d/d<)J(t) = (V •v)J{t), for the Jacobian, 
J{t) = dV^/diVo = d{x^,y*,z*)/d{xo,yo,zo), of the transformation from initial coordi¬ 
nates of fluid elements to the corresponding coordinates at time t. 

Likewise the boundary 9S7o maps onto within which area elements d5'* = nJdS* 
advect according to 


—d5‘ = (V • v)dS* - (Vv) ■ d5‘ 


(3.5) 


Representing a Lagrangian map as a dynamical flow, induced by an Eulerian velocity 
held, allows connection to be established with the modern Lie algebra approach to fluid 


dynamics, Arnold & Khesin (1998). However, although Lie operator methods are useful 


in Hamiltonian perturbation theory, Dewar (1976), the simpler Lagrangian-based varia¬ 


tional approach used in this paper avoids the need for most of this abstract machinery. 
(Likewise for abstract differential geometry.) Nevertheless we shall find the extension of 
our Eulerian-Lagrangian mapping notation to include flows other than time evolution 
makes for a compact notation, and, being semi-Eulerian, leads to a more familiar form 
for the perturbation expansion of the Lagrangian than the strictly Lagrangian approach, 
Dewar ( |1970 ). 

Application of Hamilton’s Principle requires us to vary trial fluid element pathlines 
to And the Euler-Lagrange equations that determine which such pathlines are actually 
physical. Thus we introduce a new flow that maps the position vectors of fluid elements 
from their unvaried positions at each time t to their varied positions at the same time by 


using a flow analogous to the Lagrangian map defined in eq. (3.1), but with t replaced 


by a dimensionless variation parameter e (typically small), and with the variational map 
generated by a variational velocity u. Thus, in eqs. (3.1 - |3. 5 ) replace t with e, the initial 
time tg with 0, the time evolution flow r*\v\ with the e-flow and the velocity v(x, t) 

with v>{x,e\t), to give 


dr'^{x) 

de 


= v{H{x),e) , 




(x) = X 


(3.6) 


with the parametric time-dependence (... \t) now left implicit. 

The varied Lagrangian map r*[VrJ\ is now found by composing the unvaried Lagrangian 
map with the variational map. 


r^\Vr. 


= r \u\or m 


(3.7) 


where o denotes composition of functions: f og{x) = f{g{x)) and subscript ^ denotes a 
varied quantity. This implicitly defines the varied velocity field v,^, which is to be found 
in terms of the unvaried position of a representative fluid element x^ = r*(a;g) as the 
total time derivative of the varied position, x^^ = r*[D,^](a;g) = r'^{x*\t), as 


flx^,t) = Dt a:„ 


(3.8) 


where Xr^ here denotes the Eulerian representation of the varied position, x^(x,t) = 
r’^{x\t), and Dt = dt d- v{x,t) • V is the advective derivative. 

In MRxMHD we also need to consider the case of fluid elements on the common 
interfaces (current sheets) Tij = dD,i separating subregions and Vlj. The shape 
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of the interface, which we represent as the level surface f{ij) = 0 of an appropriate 
smooth function f{x) changing monotonically across the surface, is not known a priori 
so must be subject to variation in applying Hamilton’s Principle. Thus we must introduce 
the variation parameter e in representing Tij geometrically: /(^ e) = ±0. We use 

the notation to distinguish which side of the interface a varied fluid element is on: 
/(i,j)(3^L±|^; e) = iO. Taking the total derivative of left and right sides of this expression 
with respect to t, on both sides of the interface, we find n* = 0, n = ^f{ij)/\^f(ij)\ 
being the unit normal and !•] denoting the jump in a quantity as the evaluation point 
crosses the interface (so \dtf{i^j){x = 0, as f(ij) is assumed smooth). Similarly, 
total differentiation with respect to e gives 


= 0 , 


* e r± 


(3.9) 


This states that the normal component of is constrained to be continuous across the 
interface, but otherwise it is unconstrained. Rather, T^j- is advected with i/ during vari¬ 
ations at constant t, just as it is under time evolution. 

An exception is the case where the plasma is confined by a prescribed, though possibly 
time-dependent, boundary “wall” (w). This is the special case that is not a 

function of e, so one obtains the constraint n • i/ = 0 at constant t on T)^ with the 
tangential components unconstrained. 

So far we have treated e as a finite parameter, on a par with t. However, in this 
paper we use the variational transformation only for calculating the first variation of 
the action in Hamilton’s Principle, so we need the variational map only to linear order: 
r'^\v]{x) = x + ei>{x, 0|t) -|-0(e^). Thus, defining the Lagrangian variation in position Aa; 
through x^ = x + e/S.x(x^t)+0{e^), we have Ax{x,t) = i/(a;,0|t). Similarly, dehning the 
Lagrangian variation in velocity, Av, to be such that v^{x^,t) = v(x,t) -|- eAv{x,t) + 
O(e^), we have, from eq. (3.8), 


Av{x,t) = Dt Ax{x,t) . (3.10) 

In the above we have followed [Newcornb (1962) in using A/ = limc_>o[/..^(3:~, e|t) — 
f{x^e\t)]/e to denote the Lagrangian variation in an arbitrary field /, while the corre¬ 
sponding Eulerian variation 5f is defined by 6f = lime_>o[/..^(a:, e\t) — f{x, e|t)]/e, i.e. with 
both the varied and unvaried field evaluated at the same, unvaried, position x. Thus A 
may be regarded as the operator lim£_,,o d/de while <5 is the operator limg_,,o d/de. 

By definition the two operators are related by A = 5 -|- Ax • V. Applying both sides 
to X it is easily verified, as a consistency check, that Sx = 0. Also, both S and A being 
differential operators, the product rule, e.g. A{fg) = {Af)g + gAf, and commutation 
relation 


AV/ = VA/ - (VAa;) • V/ 

apply, and correspondingly for S, where f{x,e) and g(x,e) are arbitrary. 


(3.11) 


3.2. Holonomically constrained and free variations 

In this subsection we extend the general formalism developed above for fluid kinematics 
to treat variation and perturbation of fields, specifically mass density p{x,t), pressure 
p{x,t), and magnetic vector potential A(x,t) [and hence magnetic field B = Vx A]. 

The density (mass conservation) equation is the lowest and most robust one in the 
hierarchy of moment equations used in deriving fluid models from kinetic theory. In fact 
freezing mass into fluid elements seems fundamental to any fluid theory, so we build the 
holonomic mass conservation constraints Pr..{xi,,t\e)dVf = p{x*,t)dV* = p{xo,to)dVo 


into both IMHD [as in Newcomb (1962); Dewar (1970)] and into our new MRxMHD 
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formulation. From eq. (3.4) and its e-flow analogue, these imply 


dp _ _ 


(3.12) 


or, decomposing d/dt as dt + v • V, d/de as -t- !>» • V, 

(3.13) 

As explained above, the Lagrangian variation in density is Ap = dp,..,/de|e=oi ^'iid the 


I = -V.W.|^ = -V. 


Eulerian variation is Sp = dpr^/de\,:=o. That is, from eqs. (3.12) and (3.13) 


Ap = —pV • Ax ^ 5p = —V • (pAx) . 


(3.14) 


In MRxMHD [unlike IMHD, Newcomb (1962); Dewar (1970)] pressure and magnetic 
held V X A are not holonomically constrained by the freezing of entropy and flux micro¬ 
scopically into each fluid element but are freely variable, expressed by writing pressure 
and vector potential as p{x,t\e) and A{x,t\e), with their e-derivatives 6p and SA being 
arbitrary variations in Hamilton’s Principle (as is Aa;). (This is not of course to say 
Lagrangian variations of free, i.e. not holonomically constrained, helds do not exist, but 
rather that their Eulerian variations 6 are primary, with their Lagrangian variations be¬ 
ing defined by the rule A = 5 -I- Ax • V above.) The relation of interface kinematics to 
the holonomic tangential-S constraint, eq. (1.1), is developed in Appendix [ b| 

Consider an integral of the form L = Jq CdV, where £ is a function of various helds to 
be varied in Hamilton’s Principle (in our case p, v, p and A), which are thus functions 
of e. By differentiating L with respect to e, using the e-how analogue of eq. (3.4) (after 

where, 


as m 


Frieman & Rotenberg| (1960), superscript 0 


changing variables to x°,y^, 
refers to e = 0) and integration by parts (Gauss’ theorem) we hnd the convenient identity 

6L= [ 5CdV+ [ CAx-dS, (3.15) 

JQ Jan 

which allows us to use the Eulerian variations S = dehned above to restrict the 
occurrence of Ax to those helds {p and v) that are constrained to vary with huid motions. 
[An alternative, but equivalent, approach is to transform to xq, yo, Zq at the outset, using 
Lagrangian variations and the e-how analogue of eq. (3.5) to get surface terms.] 

In applying eq. (3.15) to compute the action variation 6S = f dt f^^dVA, further 


surface terms will arise after straightforward integrations by parts to remove spatial 
derivatives of Aa;. However, the integration by parts arising from terms containing time 
derivatives is more subtle because it involves a change of variables from x, y, z to Xo,yo, zg 
to enable integration by parts, and then a change back to x, y, z. In our case the only 
such term arises from the kinetic energy term, a function of v, the time derivative arising 
from the Eulerian variation Sv = Av — Ax ■ Vv = D^Aa; — Aa; • Vv by eq. (3.10). 


Consider a Lagrangian density of the form C{v). Then, noting that Aa; is always 
taken to vanish at the endpoints of the time integration in Hamilton’s Principle, and 


using eq. (3.4) we hnd the contribution to the action from Sv, 


which has no surface term. 


dL 


dt J dEAa; • (/V • v + Vv + /D*) • — , 


(3.16) 


3.3. Lagrangian formulation of MRxMHD 

To derive the dynamics of multi-region relaxed plasmas we use the same Lagrangian 
density (1.2), integrating over each subvolume and augmenting with 
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appropriate Lagrange multiplier terms to form the Lagrangian in Li. We then sum 
to form the total effective Lagrangian 


L = L, 




(3.17) 


iGTZ 


where TZ denotes the set of plasma relaxation regions, with their Lagrangians Li being 
given by 


L,. = L 


MHD 


Tii^Si 5*20) “f Li {^i ^io'j : 


(3.18) 


where L^^^ is given by the integral of defined in eq. ( |l.2[ ), over Lti, the Si are the 

entropy invariants {mSen in the notation of Appendix]^ given by integrating the RHS 
of eq. (A4) over Lli, 


Si = 


P 


In 


K— 1 dR 


(3.19) 


da. 7 - 1 

and the magnetic helicity invariants Ki are as in eq. (2.11, evaluated over VLi. (The 
constant k, required to make the argument of In dimensionless, is arbitrary for our 
purposes but is identified physically in Appendix [A|) 

The Lagrange multipliers and ^i are constant during variation in Hamilton’s Prin¬ 
ciple (i.e are independent of e) but may change with t, as they must be chosen to make 
Si = Sio and Ki = KiQ during evolution under the Euler-Lagrange equations. The con¬ 
stant reference values S^ and Ki^ are the respective initial values at t = t^ evaluated 
over LtiQ, making Lj = when the and Hi are adjusted to satisfy the conservation 

conditions above. Subtracting off the constant KiQ also makes Li independent of the 
constant k (because of the holonomic conservation of fpdV), so the curious fact that 
the physical value of k, eq. (A5), involves Planck’s constant is irrelevant. 


If there is a vacuum region, Hv, between a plasma-vacuum interface and the wall, then 
this may be treated similarly, but with p set to zero and with the entropy and helicity 
constraints deleted, leaving the Lagrangian density 


C. = - 


BB 

2po 


(3.20) 


After applying eqs. (3.15-3.161 and appropriate integrations by parts the variation of 


the action must be of the general form 






5 A dp 
Ax , 




6x J 


(3.21) 


dQi 


where TZ'^ denotes the set of plasma relaxation regions plus the vacuum region (though 
note that only the variational derivative dS^/dA is non-zero in Hv). The interval over 
which the time integral is taken does not need to be specified as variations at the end¬ 
points are taken to vanish in Hamilton’s Principle. 

Hamilton’s Variational Principle is the statement that the Euler-Lagrange equations 
following from requiring (55 = 0 for all variations of the independent fields determine 
these fields physically. In the above variational derivative notation these are the Euler- 
Lagrange equations dSi/dA = 0, dSi/dp = 0, dSi/dx = 0 within the volumes 12^. To find 


the surface Euler-Lagrange equations the constraint eq. (3.9) needs to be considered, the 
consequences of which will be discussed in Sec. |3.5| 
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3.4. Volume variations 


Inserting eq. (3.18) in eq. (3.17), varying 5 = f Ldt, integrating by parts and comparing 
with eq. (3.211, we identify the variational derivative with respect to the vector potential 
as 

SS, 1 


(3.22) 

giving as the corresponding Euler-Lagrange equation, SSi/6A = 0, the Beltrami equation 


6A po 


eq. (2.3). In the case of the vacuum Lagrangian Lv, Vi is set to zero, giving the statement 


that the vacuum field is harmonic, V X S = 0. (While the Euler-Lagrange equation in 
the vacuum region is the same as the Beltrami equation with /r set to zero, note that we 
have actually deleted the constraint of constant magnetic helicity in the vacuum region, 
so the vacuum is not completely equivalent to a currentless plasma.) 

The variational derivative with respect to the pressure is 


Sp 


dp 


1 


7-1 


l-rA 


the corresponding Euler-Lagrange equation being 


p = np, 


(3.23) 


(3.24) 


with Ti identified [see eq. ( |A 2[ )] as the specific temperature Ti/m in 11^, where m is the 
effective ion mass (r^ is also where Ci is the ion sound speed). 

Note that, despite using internal energy and entropy densities consistent with a mi¬ 
croscopically isentropic equation of state, we have only enforced macrosopic entropy 
conservation over a whole subregion, leading to the microscopically isothermal equa¬ 


tion of state eq. (3.24) (r^ being spatially constant). Rapid equilibration of temperature 
within rij is compatible with magnetic surfaces being destroyed within a relaxation re¬ 
gion, hence poor local thermal confinement, but the bounding interfaces are assumed to 
be thermal transport barriers so the temperature can jump across each interface. Thus, 
while temperature profiles are restricted to being piece-wise constant they are otherwise 
arbitrary, enabling the use of MRxMHD to model hot, magnetically confined plasmas 


[see e.g. Hudson et al. (2012)]. 

As already remarked, the Lagrange multipliers pi and Ti may change with time to 
maintain the constancy of their respective constraints. Now that the Ti have been iden¬ 
tified as temperatures, a simple thought experiment makes it physically clear that this 
must be so: Suppose the pi, and thus the pi, are spatially constant within each then, 
just as for adiabatically deformed bags of ideal gas, TiVfi~^ = const—any change in 
volume Vi leads to a change in temperature t^. 


Using eq. (3.14) and eq. (3.16) we then find the variational derivative with respect to 


fluid element positions 


SS^ 

6x 


= -dtipv) - V • (pvv) - pV— -f pV 


9A 

dp 


Vp 


= -dtipv) - V • (pvv) -j- V ( — - 7 

7-1 \p p 

= -dtipv) - V • ipvv-i-pl) , 


(3.25) 


where the last line follows from the isothermal equation of state eq. (3.24). The corre¬ 


sponding Euler-Lagrange equation SSi/5x = 0 is the equation of motion for a compress¬ 
ible Euler fluid in momentum conservation form. Note the extraordinary simplicity of 
this result, with v being decoupled from B as the Beltrami equation eq. (2.3) implies 
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the Lorentz force j X B is zero, showing that MRxMHD supports only steady flows 
and sound waves, with phase velocity within the relaxation regions. (However 

surface waves on the interfaces can involve perturbations oi B.) 


3.5. Surface variations 

It is shown in Appendix [B| that tangential B at the interfaces Tij implies a holonomic 
constraint on variations of the tangential component, A^gt = (/ — nn) • A, of the vector 
potential, 

(M)tgt = iAxxB + V6x)tst , (3.26) 

for all X on the interface or boundary, where Sx is an arbitrary single-valued gauge poten¬ 
tial. The normal component n-SA is unconstrained. It is also shown in Appendix [B] that 
this constraint implies invariance of line integrals of A around loops on these surfaces, 
and thus conservation of flux does not need to be imposed as an extra constraint. 

Taking into account eq. (3.51 the surface variational derivative in eq. (3.21) is found to 


be 


Sx 


= \P 


ant 


\ 

- ni 

2^0 J 


(3.27) 


which is the velocity-independent part of the stress tensor, Dewar ( |1970 ) dotted with . 

However, before we can apply this result we need to take into account the fact that 
a boundary dfli is made up of interfaces between and neighbouring regions flj (say) 
across which, by eq. (3.9), the normal components of Ax are continuous. Taking into 
account the outward normals of contiguous regions being oppositely directed, this implies 
the constraint Ax-rii = —Ax-rij on the common interfaces Tij. This coupling gives the 
surface Euler-Lagrange equation, the natural boundary condition between interfaces, as 
the continuity condition 

P+S~ = ® ’ 

2/io] 

which is the same as the jump condition for advected discontinuities in ideal MHD [see e.g. 
Sec. 5.12 of Hosking fc Dewar| ( 201^ ] . (The tangential components Axtgt are separately 
variable, but give no natural boundary conditions because they do not appear in SS.) 


(3.28) 


4. Conclusion 


We have built a general framework on which to develop relaxed-MHD dynamics further. 
Some avenues to be explored are indicated below: 

• Stationary states with flow: If one invokes a modified form of the “imaginary ex¬ 
periment” of Kruskal & Kulsrud (1958) in which the fictitious friction force acts only on 


interface movements, and is sufficiently strong so as to allow only movements slow com¬ 
pared with a characteristic sound transit time, then negligible sound wave energy will be 


excited (see Sec. 2.2) and the system will relax to a static equilibrium state or one with 
steady flow. Thus action extremization would seem to provide a more physically intuitive 
framework for variational construction of equilibria with flows than one based on energy 
minimization, which requires an arbitrary angular momentum constraint, [Dennis et al.\ 
(2014) to keep the kinetic energy from being minimized to zero. 


It also does not seem necessary to invoke the cross-helicity invariant, Hameiri (2014); 


Dennis et al. (2014), mentioned in Sec. 2.1 


• Spectral and stability studies with and without flow: By including the kinetic energy 
in a natural way, our dynamical formulation of MRxMHD provides a physical normal- 
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ization for the linear growth rates of instabilities to replace the artificial one derived 


previously using an energy principle, Hole et al. (2007); Mills et al. (2009). 


The dynamical formulation also suggests performing simulations using the water-bag ap¬ 
proach, Potter (1976), for exploring the nonlinear evolution, and possible saturation, of 


instabilities. This may, for instance, help resolve the paradox that existence of equilib¬ 
rium interfaces seems to be contingent on having highly irrational rotational transforms 
on the domain boundaries dD.i, McGann et al. (2010); McGann (2013), but fixing rota¬ 


tional transform is incompatible with the constancy of the helicity invariants in general. 
Simulation may also be a useful way to explore formation of singularities at resonant 
surfaces, Loizu et al. ( 2015a|6 ). 

• Exploration of the possibility of including reconnection in MRxMHD by allowing a 
slow leak of flux and plasma through the interfaces. 


Appendix A. Application of ideal-gas thermodynamics to plasmas 

The thermodynamics involved in MRxMHD is elementary, being the same as for an 
ideal gas. However, the expressions used here for the internal energy density p/ij — 1) 


and entropy constraint density plii{Kp/p'^)/— 1) introduced in Bhattacharjee & Dewar 
(1982) (an essentially arbitrary constant quantity k here being inserted to make the argu¬ 


ment of the logarithm dimensionless) are somewhat different from the expressions found 
in most thermodynamics texts. Thus we briefly review their derivation from standard 


thermodynamics and its adaptation to MHD (extending the discussion in Dewar et al. 

( [2008t ). 

First recall that, for a single-species ideal gas of absolute temperature Tk (in degrees 
Kelvin) whose atoms are of mass m and number density is n, the mass density p is mn 
and the pressure p is nk-aT]^, where fee is Boltzmann’s constant. The internal energy U 
is {3/2)nVkBTK = pVK^ — 1), where V is the volume of the system and 7 = 5/3 is the 
ratio of specific heats. The statistical mechanical entropy Sk (in units such that a heat 
increment is dQ = TxdS'K) is given by the Sackur-Tetrode equation 


S'k = IV/cb < In 


V Mirm 
N \^n) 


(Al) 


where N = nV is the number of particles and h is Planck’s constant. 

In plasma physics, temperature T is measured in energy units, i.e. T = k^T^, the 
corresponding entropy in energy units being ^en = in order that dQ = Td^en- 

Also there are two species, ions and electrons, to take into account, their number densities 
being denoted n-, and rig, respectively. If is the effective ionization state then, to a very 
good approximation, Uj = n^/Zef^ to maintain quasineutrality. Then the total pressure 
p = UgTg + n{Ti becomes n^iT^ Ti/Zgff). 

The MRxMHD assumption that current sheets on magnetic surfaces act as transport 
barriers is most justifiable if T,, the small gyroradius of the electrons providing good 
confinement across magnetic field lines and their rapid motion along field lines providing 
fast thermal equilibration on magnetic surfaces and within the chaotic relaxation regions. 
Thus we henceforth assume the ions are cold, Ti/Tg <C 0. However, the mass density is 
dominated by the ions because TOg/wi ^ 0 . 

Simplifying notation by denoting ng by n and Tg by T, and defining an effective particle 
mass m = milZf,Q, we summarize these approximations as 

(A 2) 


p = nm and p = nT . 
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To adapt standard thermodynamics we model the plasma as a monatomic gas at 
temperature T made up of particles of mass m [except in the de Broglie term in eq. (A 1), 
where we use rrie] and write the Sackur-Tetrode equation eq. (A 1) as 

mSen = VSm , (A 3) 


where we have derived the entropy constraint density used in eq. (3.18) as 

P 


7-1 


■ In ( 


the hitherto arbitrary non-dimensionalizing constant k now being identified as 

AirmeimeV 
3(7-l)h2 ■ 


(A 4) 


(A 5) 


Appendix B. Vector Magnetic Potential Boundary Constraints 


In this appendix we seek to justify the holonomic constraint eq. (3.26) on a plasma or 


vacuum region boundary dVL. Also, to verify that magnetic fluxes are conserved under 
variation, we need to show line integrals j>A-dl around loops on the interface are invariant 
under displacements of the interface. While eq. ( 3.26|) is as expected from the ideal MHD 


result, Bernstein et al. (1958 1, that SB = Q = V X {Ax X B), it needs to be justified for 


MRxMHD because we make no frozen-in-flux assumption other than the tangential-B 
constraint. Within the subregions the IMHD result SB = Q does not in general 
apply. 


We first consider the problem of propagation of the tangential-B condition eq. (1.1) on 
a time-dependent surface T* and then adapt the results to find the analogous e-variations 


at fixed t. First, from eq. (1.1) and eq. (|3.5|) we have (suppressing the superscripts t unless 


needed to emphasize time dependence) 


dS-B-{Vv)-dS = 0 


(Bl) 


(B2) 


Dividing by dS and using the definition dB/dt = dB/dt + v-VB we thus find 

dB 

n ■ = n • V X (u X B) on T . 

(NB This is obviously consistent with the IMHD equation dtB = Vx(uxB), but is 
derived completely generally and is thus applicable to MRxMHD as well.) 

Substituting B = V X A in eq. (B 2 1 we easily find n • V X {dtA — u X B) = 0, which 
is equivalent to 


|V/| 


V/x 


dA 

Ik 


-vxB 


= 0 on T 


(B3) 


where / is a differentiable function such that f{x,t) = const on T* (as in the last two 


paragraphs of Sec. 3.1) and we have used the identity V X V/ = 0. 


Although eq. (B3) is written using 3-dimensional Cartesian vector calculus notation. 


it applies only on the 2-dimensional surface T. We now resolve this seeming paradox 
by transforming to a curvilinear coordinate system /, 5 , h such that the basis vectors 
= V/ = n|V/|, e® = Vg, and = V/i are linearly independent, so that the pair 
h, g specifies a point on T: / = const and V = e-^df + e^dg + e^dh- 
Using the identity >7 V • u = df{Jef‘u)-Pdg{Je^‘u)-\- dh{Je^ • u), where J = 
l/g/ . e® X and u is an arbitrary vector field, and choosing u = V/ X {dtA — vx B), 
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we find eq. (B 3) is equivalent to 


dn 


JVg-Vfx 


dA 

iR 


-vxB 


+ dh 


JVh ■ V/ X 


dA 

iR 


-vxB 


= 0 . 


(B4) 


As the left-hand side does not contain the normal derivative df, it is a surface divergence 
operating purely on values of u evaluated at the surface F. Thus the nature of the spatial 
dependence of B{x,t) off F* is immaterial to the evaluation of the boundary condition 
eq. (B3)—in particular B does not need to have nested magnetic surfaces. [Likewise 
the off-surface dependence of f{x, t) is irrelevant, as | V/| cancels in the product JVf = 
n/n-VgxVh.] We also note, using eq. (1.11, that dtA — vxB = Vfx[{dtA)tgt+VnB], so 
that only the tangential components of dtA and the normal velocity component = n-v 
contribute. 


Clearly, the general solution of eq. (B 3) is V/X (9*A — vxB) = V/ xVdtXi where x 
is an arbitrary gauge potential. Crossing both sides with n/jV/j and rearranging gives 
the alternative form 

. (B5) 

^^Agt V ^Atgt 

Replacing t with e and v = dr* (x)/dt with dr*^ (a:)/de in eq. (B 51 and taking the limit 
as e —>■ 0 [cf. discussion after eq. (3.10)] gives the desired variational holonomic constraint 
eq. ( 3.26| ) . 


To show invariance of loop integrals j>A-dl under boundary and interface displace¬ 
ments, on surfaces that are not simply connected, we first show A(A • dl) — AA ■ di -|- 
A ■ AdZ is a complete differential on these surfaces. From the epsilon-flow analogue of 
eq. (3.3), AdZ = dZ • VAa;. Using the constraint eq. (3.26) we find 


A{A . dZ) = dZ • [6A + AX-VA+ (VAa;) • A] 
= dl • [Aa; x (V X A) -f VSx 

-l-Aa; • VA -I- (V Ax) • A] 

= dl ■ V{Ax ■ A -k V6x) , 


(B6) 


which is a perfect differential as required. Thus, there is zero variation in line integrals 
around loops provided we also require Sx to be single-valued. 
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